Finite-momentum Bose-Einstein condensates in shaken 2D square optical lattices 



M. Di Liberto 1 ' 2 , 0. Tieleman 1 , V. Branchina 3 ' 4 , and C. Morais Smith 1 
1 Institute for Theoretical Physics, Utrecht University, Leuvenlaan 4, 3584CE Utrecht, the Netherlands 
2 Scuola Superiore di Catania, Universitd di Catania, Via Valdisavoia 9, 1-95123 Catania, Italy 
3 Department of Physics, University of Catania and 
A INFN, Sezione di Catania, Via Santa Sofia 64, 1-95123 Catania, Italy 
(Dated: January 20, 2013) 

We consider ultracold bosons in a 2D square optical lattice described by the Bose-Hubbard model. 
In addition, an external time-dependent sinusoidal force is applied to the system, which shakes the 
lattice along one of the diagonals. The effect of the shaking is to renormalize the nearest-neighbor 
hopping coefficients, which can be arbitrarily reduced, can vanish, or can even change sign, depending 
on the shaking parameter. It is therefore necessary to account for higher-order hopping terms, which 
are renormalized differently by the shaking, and introduce anisotropy into the problem. We show 
that the competition between these different hopping terms leads to finite- momentum condensates, 
with a momentum that may be tuned via the strength of the shaking. We calculate the boundaries 
between the Mott-insulator and the different superfluid phases, and present the time-of- flight images 
expected to be observed experimentally. Our results open up new possibilities for the realization of 
bosonic analogs of the FFLO phase describing inhomogeneous superconductivity. 

PACS numbers: 03.75.-b, 03.75.Lm, 67.85.-d, 67.85.Hj 



I. INTRODUCTION 

Ultracold atoms in optical lattices are ideal systems 
to simulate the quantum behavior of condensed matter 
because the lattice geometry, the type of atoms (bosons 
or fermions), and their interactions can be manipulated 
in a perfectly clean environment pQ. Furthermore, they 
provide a perfect testing ground for a wide variety of 
theoretical models. One of the most prominent exam- 
ples is the Bose-Hubbard model, which has been studied 
extensively theoretically (e.g. Refs. [2j [3]) and realised 
experimentally [4] [5] . 

More recently, much interest has been devoted to time- 
dependent, periodically stirred optical lattices, which al- 
low for engineering synthetic gauge fields into the system 
[6HH]. In the presence of a staggered rotation, Dirac cones 
were shown to emerge for square optical lattices, thus 
simulating the behavior of anisotropic graphene when the 
system is loaded with fermions [9 . If instead a mixture 
of fermions and bosons is used, several properties of high- 
T c superconductors can be reproduced [TO] [11] . Loading 
the same lattice with dipolar bosons leads to a supersolid 
phase with vortices [T2] . 

Besides the interesting features that arise in the pres- 
ence of rotation, a full range of new possibilities was 
shown to emerge by shaking the optical lattice. If the 
shaking frequency is much larger than the other char- 
acteristic energy scales in the problem, the parameters 
of the Hamiltonian are renormalized. This provides an- 
other tool to control the lattice parameters, and even en- 
ables the simulation of otherwise experimentally inacces- 
sible lattice models [T3HT5] . In the Bose-Hubbard model, 
for instance, the superfluid-Mott-insulator transition has 
been driven by ramping the shaking perturbation and 
thus tuning the effective hopping parameter to zero [16] . 
Another fascinating experiment has revealed that mag- 



netically frustrated systems can be realized with spinless 
bosons by applying elliptical shaking to a triangular lat- 
tice [17]. 

Here, we consider a 2D square lattice shaken along 
one diagonal and investigate the effect of next-nearest- 
neighbor (nnn) and next-next-nearest neighbor (nnnn) 
hopping in the behavior of a bosonic system. In an 
effective description, the shaking perturbation leads to 
a renormalization of the nearest-neighbor (nn) hopping 
parameter, which can vanish or even become negative 
[13] [15]. When this parameter is tuned to be very small, 
higher order hopping terms, which are usually negligible, 
may become relevant and must therefore be included in 
the model. Although the nnn hopping coefficients are 
strictly zero in 2D optical lattices where the x— and y— 
directions are independent (separable potential) , they are 
relevant for non-separable optical lattices. In this paper, 
we show that a tunable finite-momentum condensate can 
be realized in a certain range of parameters for a realistic 
and simple setup, thus bringing us a step further in the 
realization and control of finite- momentum Bose-Einstein 
condensates (BECs). 

Finite-momentum condensates have recently attracted 
a great deal of attention. In the original proposals by 
Fulde, Ferrel, Larkin, and Ovchinnikov (FFLO), it was 
argued that finite momentum Cooper pairs would lead 
to inhomogeneous superconductivity, with the supercon- 
ducting order parameter varying spatially (the so-called 
FFLO phase) [18 . Early NMR experiments at high mag- 
netic fields and low temperatures in the heavy-fermion 
compound CeCoIns have shown indications of an FFLO 
phase [19], although recent data suggest the existence of 
a more complex phase, where the exotic FFLO supercon- 
ductivity coexists with an incommensurate spin-density 
wave [20] . For ultracold fermions with spin imbalance, 
on the other hand, the observation of the FFLO phase 
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has been recently reported in ID [2T] . 

Earlier theoretical studies of a square-lattice toy model 
for a scalar field, which took into account non-trivial hop- 
ping beyond nearest neighbors, have shown that quantum 
phases may be generated in which the order parameter is 
modulated in space [22]. Finite-momentum condensates 
were also experimentally detected for bosons in more 
complex lattice geometries, such as the triangular lattice 
under elliptical shaking [17] , or for more complex interac- 
tions, as e.g. for spinor bosons in a trap in the presence 
of Zeeman and spin-orbit interactions [23 . With regard 
to bosons in a square lattice, it was recently shown that a 
staggered gauge field may lead to finite-momentum con- 
densates [24 . In this case, the bosons condense either at 
zero momentum or in the corner of the Brillouin zone, 
and a first-order phase transition occurs between these 
two phases [24]. Here, we propose that finite- momentum 
condensates can be realized for bosons in a shaken square 
lattice and that we may tune the momentum of the con- 
densate smoothly from to (7r,7r), by varying the shak- 
ing parameter Kq. To the best of our knowledge, this is 
the first time that such an effect has been predicted for 
optical lattices as originating solely from beyond-nearest- 
hopping terms. The interaction simply shifts the ground- 
state energy by a constant and does not change the con- 
densation momentum. 

In the following, we consider in Sec. II an ex- 
tended Bose-Hubbard model which includes higher- 
order-hopping coefficients for a non-separable 2D square 
optical lattice. After introducing a sinusoidal shaking 
force to the system, we show in Sec. Ill how the finite- 
momentum condensate arises, and how the condensation 
momentum depends on the shaking. We present a 3D 
phase diagram, with as parameters the Hubbard inter- 
action J7, the chemical potential /i, and the shaking pa- 
rameter i^o, and indicate the required parameters for the 
realization of the tunable regime in Sec. IV. Finally, we 
calculate the expected outcome of time-of-flight experi- 
ments in Sec. V and present our conclusions in Sec. VI. 



II. THE MODEL 

Before discussing the generic 2D problem, let us recall 
the behavior of ID lattices and 2D separable lattices. A 
simple calculation shows that in ID optical lattices of the 
form V(x) = (Vb/2) cos(2fcx) (Vb is the potential depth 
and k = 2ir/\ is the wave vector of the laser beam), 
nnn hopping coefficients do not change the position of 
the global minima in the single-particle spectrum but 
generate metastable states. 2D separable potentials do 
not introduce new physics from this point of view. The 
simplest non-separable potential in 2D is given by [25 

V (r) = -V { sin 2 [k(x + y)\ + sin 2 [k(x - y)\ 

+ 2asin[k(x + y)] sm[k(x — y)] f , 



where k = 2v / 2tt/2 A (the factor a/2/2 comes from a coor- 
dinates transformation corresponding to a rotation of the 
lattice of 7r/4) and we will make the choice a = 1 in the 
remainder of this work. Had we chosen a = 0, the poten- 
tial would have been separable, whereas for < a < 1 
the potential would correspond to a superlattice, with 
neighboring wells of different depths. 




FIG. 1. Non-separable optical potential V(x,y) given by Eq. 
|l]) with a = 1. 

As shown in e.g. Ref. pQ, we can calculate the hopping 
coefficients from the exact band structure 

£„(q) = ^ n (R)e^ R , (2) 

R 

where n is the band index, q is the quasimomentum, and 
R is a lattice vector. In this notation, t n (H) is the hop- 
ping coefficient between two sites separated by the lattice 
vector R in the n-th energy band. The non-separable op- 
tical potential generates hopping coefficients along direc- 
tions other than those of the elementary lattice vectors of 
the lattice which were exactly zero for separable poten- 
tials. A lattice vector has the form R = mae x + nae yi 
where a = X/V2 is the lattice spacing, m and n are 
integers, and e x and e y are unit vectors in the x- and y- 
directions; R is indicated in short notation as (m, n). For 
the non-separable potential that we have introduced, we 
find non zero hopping terms also for pairs of sites identi- 
fied by (1,1) or (2,1), which vanish for separable lattices. 
Table[l] shows the most relevant lowest-band hopping co- 
efficients for shallow lattices. Higher-order hopping coef- 
ficients are neglected because they are at least ten times 
smaller than t" and therefore not important, as will be- 
come clear afterwards. 



Vo/E r 


(1,0) -t 


(1,1)***' 


(2,0) t" 


1.0 


-2.45 x 10" 2 


-8.89 x 10" 4 


8.88 x 10" 4 


2.0 


-4.52 x 10" 3 


-6.65 x 10~ 5 


2.27 x 10" 5 


3.0 


-1.06 x 10" 3 


-5.89 x 10~ 6 


1.06 x 10" 6 


4.0 


-2.97 x 10" 4 


-6.74 x 10" 7 


7.86 x 10" 8 



TABLE I. Relevant hopping matrix elements (in units of the 
recoil energy E r ) of the lowest band for shallow lattices. 
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We will assume that the lowest-orbital Wannier func- 
tions are still even and real for this non-separable poten- 
tial. As shown by Kohn [26], this can be proven for sep- 
arable potentials; for non-separable ones it is also a rea- 
sonable conjecture, supported by numerical simulations, 
as shown in Ref. [27]. If we apply a driving sinusoidal 
force like the one studied in Ref. [13] , but now along one 
of the diagonals, the shaking term in the co-moving ref- 
erence frame that has to be added to the Hamiltonian 
reads 

W(t) = K cos(u;t) ^(i + j)n^ , (3) 



where uj is the shaking frequency, r is the real time, and 
riij is the density operator at site (i, j). Following the 
approach discussed in Refs. p~3] [15], the non-interacting 
effective hamiltonian for the quasienergy spectrum in the 
high-frequency limit hu ^> U,t (and thus hu ^> t',t") is 



H ef£ = -tJo(Ko) ^ a i a r±*v +^0(2if )^4 a r±(e x +e y ) +i'^4 fl r±( ei -e„) + t" Jq(2Kq) ^ 4 a r±2e^(4) 
r,v—x,y r r r,u—x,y 



where the shaking parameter is Kq = K/hw. The Bessel 
function Jo(x) has a node at x ~ 2.4048; hence, when 
the nn- hopping coefficient t e ff = t Jo(Ko) is negligible, 
the higher-order ones are not. Note that the hopping co- 
efficient along the diagonal perpendicular to the shaking 
direction is not affected by the shaking. 



III. TUNABLE FINITE-MOMENTUM 
CONDENSATE 

The effective Hamiltonian is diagonal in reciprocal 
space and the single-particle spectrum reads 

E h = -2tJ (K ) [cos(k x ) + cos(ky)} 

-h2£' J (2iTo) cos(k x + k y ) + 2t' cos(k x - k y ) 
+2t"J (2K ) [cos(2k x ) + cos(2^)] , (5) 

where k v — k • e u and we have set the lattice constant 
to unity. The spectrum has an absolute minimum at the 
center of the Brillouin zone (k = 0) when K$ < 2.4048 — 6 
and at the four corners of the Brillouin zone when Kq > 
2.4048 + 5. In the interval 2.4048 - S < K < 2.4048 + 5, 
two symmetric minima develop along one diagonal of the 
Brillouin zone at ±ko, as shown in Fig. [2] We may deter- 
mine the precise position of these two minima by study- 
ing the first derivative of the single-particle spectrum for 
k == k x — ky . The non-trivial minima are given by the 
solution of the equation 



cos(ka) 



Jo(K ) 



2J {2K ){t l +2t 2 ) 



f(K ) 



(6) 



where t\ = t 1 jt and = t" jt. We have found that for 
Vo = 2E r , 3E r , and AE r , the second derivative of the 
single-particle spectrum shows that Eq.([6| corresponds 
to a true minimum, while for Vo = lE r it is a maxi- 
mum. The largest interval E = 26 of the shaking pa- 
rameter Kq for which the non-trivial minima appear has 



Ejt/2t 




FIG. 2. Single-particle spectrum at Vo 
and contour plot. 



3E r for K = 2.4048 



been found to be at lattice depth Vo = 2.2E ri where 
5 = 0.0045, and hence the condensation momentum is 
finite for 2.4003 < K < 2.4093. Since we expect the 
bosons to condense at the minimum of the single-particle 
spectrum, the condensation momentum given by Eq. ^ 
is a function of the shaking parameter Kq and smoothly 
evolves from k = at the left edge of E to k = (7r, tt) at 
the right edge of E, see Fig. [3] The two minima in the 
E region are inequivalent because they are not connected 
by reciprocal lattice vectors and we thus need to take 
both into account for evaluating the condensation mo- 
mentum. The arccosine shape of the evolution curve can 
be explained by linearising Eq. |6| around Kq = 2.4048, 
which is a good approximation because 6 <C 1. 



The non-interacting ground state of the tunable- 
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momentum SF phase with momenta ±ko is 

N 



(<) n (al ko ) N - n \o) 



% y/n\(N - n)\ 



N 



= ^c„|n ko ,(iV-n)_ ko ), 



(7) 



n=0 



where the coefficients c n obey the normalization condi- 
tion \ c n\ 2 = 1- The ground state is thus N + 1-fold 
degenerate, where N is the number of particles. 




2.4 2.402 2.404 2.406 2.408 2.41 



FIG. 3. Evolution of the minimum in the single-particle spec- 
trum in units of the lattice spacing a as a function of the 
shaking parameter Ko at Vb = 3E r (only the positive branch 
is considered). 

We stress that there is a close similarity between our 
system and a BEC of magnons (or triplons). In dimerized 
antiferromagnets, the magnons condense at a non-zero 
wavevector ko = (tt/cl, it /a) for applied magnetic fields H 
which lie between two critical values, for H c i < H < H C 2 
[28] . In addition, the magnon dispersion of a two-leg 
antiferromagnetic ladder with frustrated nnn couplings 
along the legs shows a minimum that is incommensurate 
with the lattice spacing [29 . 



where N s is the number of lattice sites. The matrix ele- 
ment is in diagonal form and the eigenvalues are an up- 
side down parabola in n. This means that the minima 
are at the edge of the interval n £ [0, N] and that they 
are degenerate. Interactions have partially removed the 
degeneracy; the perturbative (degenerate) ground state 
to zeroth order is 

c+ p(«L) JV |o) + ^(« t - kn ) N |o) (io) 



\G) 

and has energy 



>N\ 



U 



(H) = (H ) + (H int ) = NE k0 + ^N(N - 1). (11) 

Eq. ( p~Q| ) shows that the ground state is a superposition 
of two degenerate states in which all the particles have 
momentum ko or — ko. These two states are entangled 
and behave in a very similar way to the states found 
by Stanescu et al [30] for condensates with spin-orbit 
coupling. 

One can generalize the approach described in Ref. [31] 
to calculate the MI-SF phase boundaries, taking into ac- 
count higher-order hopping terms. The outcome is 



p± = -(2N 



1)- 



2(27V + l)^ ko +c7 2 , (12) 



where ft = /x/2t, U = U/2t, £k = E\^ Q /2t and k is the 
condensation momentum, which depends on the shaking 
parameter Kq. Plotting jl± then gives the phase diagram, 
which is shown in Figs. [4] and [5] We note that the con- 
densation momentum is not changed by the interactions. 
This can be seen e.g. by doing first order perturbation 
theory calculations: in the presence of interactions, the 
energy per particle is shifted by an amount of NU/2N S , 
which is momentum independent. 



V. EXPERIMENTAL CONSIDERATIONS 



IV. PHASE DIAGRAM 

Let us now consider an additional term to the Hamil- 
tonian Q, which describes the local interactions between 
the bosons 

#int = ^^rK-1). (8) 
r 

We will treat the interactions between the atoms in a 
perturbative way and study their effect on the ground 
state degeneracy. By applying first-order perturbation 
theory, we find that the correction to the ground-state 
energy NE\^ is given by 

(m, N — m\H{ n t\n, N — n) = 

= ^-[-2n 2 + 2nN + N(N-l)]5 mn , (9) 



The lobe with unit filling No = 1 yields a critical value 
of U below which only the SF phase is allowed, see Fig. [6] 
Typical values of U are too large to allow us to probe 
the tunable-momentum SF with ordinary experimental 
setups. However, we can decrease U by reducing the s- 
wave scattering length with Feshbach resonances, which 
are available for both the Rubidium isotopes and also for 
other alkali atoms. We remark that although the range 
of Kq in which the condensation momentum is tunable 
is very small, the required precision is well within exper- 
imental control [32] . 

The quantum phases discussed above could be exper- 
imentally observed by doing the usual time-of- flight ex- 
periments. These experiments measure the momentum- 
space density distribution 

n(k) = (V f (k)^(k)) 

= N\W(k)\ 2 (|c+| 2 o k , ko + |c_| 2 o k ,_ ko ) , (13) 
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(a) 




FIG. 4. Phase boundaries for Vb = 3E r where ji = /n/2t and 
U = U/2t: (a) lobe with N = 1; (b) lobes with JV = 1, 2 in 
the region of the tunable finite- momentum condensate. 
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FIG. 5. Phase boundaries between the Mott-insulator and the 
superfluid phase for Vb = 3E r at fixed Ko and filling factor 
iV = 1: ( a ) Ko = 1; (b) K = 2.405; (c) K = 3. Note the 
different scales for each plot. 



where W(k) is the Fourier transform of the Wannier func- 
tion and we adopted the coherent state approximation for 
the SF ground state. The delta functions select the po- 
sitions of the peaks in the absorption image and are a 
clear signal of the presence of such a condensate. When 
the image is recorded, and the first atom is measured to 
have one of the two momenta, then the wave function 
collapses on that state, showing only one peak. An array 
of identical 2D systems would reveal a pattern with both 
peaks and we can study the effect of the non separability 
of the optical potential on the Wannier functions in recip- 
rocal space. In Fig|7| we show a qualitative indication of 
the time-of-flight image described by Eq.(13). We have 
lumped together the effects of a hypothetical external 
trap and the Fourier transform of the Wannier function 
into a Gaussian filter, suppressing the peak heights in 
higher Brillouin zones. In addition, we have modeled the 
broadening of the peaks by replacing every peak by a 
highly localised Gaussian. 

It is instructive to compare this pattern with predic- 
tions for other systems like the time-of-flight images given 
in e.g. Refs. [9[ [12] for finite-momentum superfluids and 
supersolids. We clearly see that these phases have a 




FIG. 6. Critical value U c = (U/2t) c at Vb = 3E r as a function 
of the shaking parameter Ko; in the inset, we show the region 
2.4016 < Ko < 2.4081 where the tunable finite-momentum 
condensate is generated. 



pattern different from Fig(7| because the position of the 
peaks differs. Hence, we can be sure that we would have 
an unambigous signal of the measurement of the tunable- 
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FIG. 7. Time-of- flight picture expected from experiment as a 
signal of the finite-momentum-condensate phase. The black 
square indicates the first Brillouin zone, and the condensation 
momentum represented in this image is ko = 27r/5, 2tt/5. 

momentum superfluid phase from experiments. 

VI. DISCUSSIONS AND CONCLUSIONS 

In conclusion, we have explored the possibility to gen- 
erate finite- momentum condensates in optical lattices un- 
der shaking, where the suppression of hopping can be 
tuned by the shaking. This opens up the possibility to 
investigate the role of higher-order hopping. To look for 
nontrivial condensation points lying inside the first Bril- 
louin zone, we have studied non-separable optical poten- 
tials in 2D square lattices. By applying the shaking along 



the diagonal of the lattice, we found that in a small region 
of the shaking parameter, where the nn tunneling is sup- 
pressed, two kinds of higher-order-hopping coefficients 
govern the dynamics of the condensate. In this region, 
we have unveiled an intermediate phase, for which the 
condensation point varies continuously from the center 
to the edge of the Brillouin zone as we tune the shaking 
parameter. There are two minima in the single-particle 
spectrum and they are symmetric with respect to the 
center of the Brillouin zone. In addition, we found that 
small interactions between the particles force the ground 
state to be a superposition of two possible Fock states: 
one where all the particles condense in one minimum and 
the other where all the particles condense in the second 
minimum. 

Finally, we note that the tunable-momentum conden- 
sate can be measured experimentally if the on-site inter- 
action is reduced significantly. This can be achieved with 
present state-of-the-art experimental techniques, and we 
hope that our results can stimulate further experiments 
in this direction. This work yields new insights for the 
realization and control of finite-momentum condensates 
and opens new possibilities in the search for bosonic 
analogs of the FFLO superconducting phase. 
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